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Abstract 

Making use of the operator product expansion, we derive a general class of sum rules for the imaginary 
part of the single-particle self-energy of the unitary Fermi gas. The sum rules are analyzed numerically with 
the help of the maximum entropy method, which allows us to extract the single-particle spectral density as a 
function of both energy and momentum. These spectral densities contain basic information on the properties 
of the unitary Fermi gas, such as the dispersion relation and the superfluid pairing gap, for which we obtain 
reasonable agreement with the available results based on quantum Monte-Carlo simulations. 
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I. INTRODUCTION 


The unitary Fermi gas, eonsisting of non-relativistie fermionie partieles of two speeies with 
equal mass, has been studied intensively during the last deeade [1-3]. The growing interest in this 
system was prompted espeeially by the ability of tuning the interaetion between different fermionie 
speeies in ultraeold atomie gases through a Feshbaeh resonanee by varying an external magnetie 
field. This teehnique allows one to bring the two-body seattering length of the two speeies to 
infinity and therefore makes it possible to study the unitary Fermi gas experimentally. Using 
photoemission speetroseopy, the measurement of the elementary exeitations of ultraeold atomie 
gases has in reeent years beeome a realistie possibility [4, 5]. Understanding these elementary 
exeitations from a theoretieal point of view is henee important and a number of studies devoted to 
this topie have already been earned out [6-8]. We will in this work propose a new and independent 
method for eomputing the single-partiele speetral density of the unitary Fermi gas, whieh makes 
use of the operator produet expansion (OPE). 

The OPE, whieh was originally proposed in the late sixties independently by Wilson, Kadanoff 
and Polyakov [9-11], has proven to be a powerful tool for analyzing proeesses related to QCD 
(Quantum Chromo Dynamies), for whieh simple perturbation theory fails in most oases. The rea¬ 
son for this is the ability of the OPE to ineorporate non-perturbative effeets into the analysis as 
expeotation values of a series of operators, whieh are ordered acoording to their sealing dimen¬ 
sions. Perturbative effeets ean on the other hand be treated as ooeffioients of these operators (the 
“Wilson-ooeffioients”). The OPE has speeifieally been used to study deep inelastio seattering pro- 
oesses [12] and has espeeially played a key role in the formulation of the so-oalled QCD sum rules 
[13, 14]. 

In recent years, it was noted that the OPE can also be applied to strongly coupled non- 
relativistie systems such as the unitary Eermi gas [15-27]. Initially, the OPE was used to rederive 
some of the Tan-relations [28-30] in a natural way [15] and, for instance, to study the dynamic 
structure factor of unitary fermions in the large energy and momentum limit [20, 23]. Eurther- 
more, the OPE for the single-particle Green’s function of the unitary Eermi gas was computed 
by one of the present authors [25] up to operators of momentum dimension 5, from which the 
single-particle dispersion relation was extracted. As the OPE is an expansion at small distances 
and times (or large momenta and energies), the result of such an analysis can be expected to 
give the correct behavior in the large momentum limit and is bound to become invalid at small 
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momenta. The analysis of [25] confirmed this, but in addition somewhat surprisingly showed that 
the OPE is valid for momenta as small as the Fermi momentum /cp, where the OPE still shows 
good agreement with the results obtained from quantum Monte-Carlo simulations [7]. 

The purpose of this paper is to extend this analysis to smaller momenta, by making use of the 
techniques of QCD sum rules, which have traditionally been employed to study hadronic spectra 
from the OPE applied to Green’s functions in QCD. Our general strategy goes as follows: 

• Step 1: Construct OPE 

At first, we need to obtain the OPE for the single-particle Green’s function Q^ik^.k) in 
the unitary limit, which can be rewritten as an expansion of the single-particle self-energy 
E|(ko,fc). The subscript '\ here represents the spin-up fermions. The main work of this 
step has already been carried out in [25]. L^{kQ,k) can be considered to be an analytic 
function on the complex plane of the energy variable ko, with the exception of possible cuts 
and poles on the real axis. Considering the OPE at T = 0, with equal densities for both 
fermionic species {n^ = n^) and taking into account operators up to momentum dimension 
5, the only parameters appearing in the OPE are the Bertsch parameter and the contact 
density, which are by now well known from both experimental measurements [31-33] and 
theoretical quantum Monte-Carlo calculations [34, 35]. 

• Step 2: Derive sum rules 

From the fact that the OPE is valid at large |ko| and the analytic properties of the self-energy, 
a general class of sum rules for ImE'[-(( 0 , k) can be derived. In contrast to the complex ko, 
CO here is a real parameter. These sum rules are relations between certain weighted integrals 
of ImE|(( 0 , k) and corresponding analytical expressions that can be obtained from the OPE 
result (for details see Section II): 

/ oo 

Jm/C(£0,M)ImE|(m,fc). (I) 

-oo 

The kernel K.{co,M) here must be an analytic function that is real on the real axis of co 
and falls off to zero quickly enough at (0 —)■ -l-oo, while M is some general parameter that 
characterizes the form of the kernel. In the practical calculations of this paper, we will use 
the so-called Borel kernels of the form 

• Step 3: Extract ImE^((0,fc) via MEM and obtain ReE^((0,fc) from the Kramers-Krdnig 
relation 
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As a next step, we use the maximum entropy method (MEM) to extraet the most probable 
form of ImE'[-((0,fc) from the sum rules, following an approaeh proposed in [36] for the 
QCD sum rule ease. 


It should be mentioned here that this method is somewhat different from the analysis pro¬ 
cedure most commonly employed in QCD sum rule studies, where the spectral function 
(which corresponds to ImE-j- here) is parametrized using a simple functional ansatz with a 
small number of parameters which are then fitted to the sum rules. This method has tradi¬ 
tionally worked well if some sort of prior knowledge on the spectral function is available and 
assumptions on its form can thus be justified. On the other hand, in cases where one does 
not really know what specific form the spectral function can be expected to have, sum rule 
analyses based on (potentially incorrect) assumptions on the spectral shape always involve 
the danger of giving ambiguous and even misleading results. MEM is therefore our method 
of choice, as it allows us to analyze the sum rules without making any strong assumption 
on the functional form of the spectral function and hence makes it possible to pick the most 
probable spectral shape among an infinitely large number of choices. 


Once ImE^(( 0 , k) is obtained from the MEM analysis of the sum rules, it is a simple matter 
to compute ReE|((0, k) by the Kramers-Kronig relation. 


/ I ,ImEt((o',fc) 

ReEt to, k) = -P / d(o' - A 

' n J-oo CO - to' 


( 2 ) 


• Step 4: Compute single-particle spectral density 

Erom the real and imaginary parts of the self-energy, the single-particle spectral density can 
then be obtained as, 

= Iro-—--———, (3) 

' n to-fzO+-efe-E|(co-l-iO+,fc) 

where is defined as = k'^/ (2m), with m being the fermion mass. 

The above steps are shown once more in pictorial form in Eig. 1. 

As a result of the above procedure, we find a two-peak structure in the imaginary part of the 
self-energy, the two peaks moving from the origin (to = 0) to positive and negative directions of 
the energy with increasing momentum \k\. Translated to the single-particle spectral density, this 
leads to a typical superfiuid BCS-Bogoliubov-like dispersion relation with both hole and particle 
branches and a nonzero gap value. 
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region of bad OPE 
convergence 


•Step 1; 

Construct OPE for the 
self-energy E(A:o,fc), which 
is valid at large |A:o| 



•Step 2: 

Derive sum rules for 
ImE(a), fe) on the real axis 
from E°™(A:o, k) at large |A:o| 


•Step 3; 

Extract ImE(a),fe) via MEM 
and obtain ReE(a),fe) by the 
Kramers-Kronig relation 


•Step 4; 

Compute the single-particle 
spectral density A(cu, fc) from 
T,{(o,k) 




FIG. 1. Steps for computing the single-particle spectral density from the OPE of the single-particle Green’s 
function of a fermionic operator. 


The paper is organized as follows. In Seetion II, we diseuss the OPE of the single-partiele 
Green’s funetion and explain how it ean be rewritten as an expansion of the single-partiele self¬ 
energy. Next, we outline the derivation of the sum rules from the OPE. In Seetion III, the MEM 
analysis results of the sum rules are shown and the consequent final form of the single-particle 
spectral density and the dispersion relation are presented. The spectral density is visualized in Eig. 
5 as a density plot and the detailed numerical properties of the dispersion relation are described in 
Table II. Section IV is devoted to the summary and conclusions of the paper. Eor the interested 
reader, we provide in the appendices detailed accounts of the relevant calculations, which were 
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needed for this work. 


II. FORMALISM 

A. The operator product expansion 

The operator produet expansion (OPE) is based on the observation that a general produet of 
non-loeal operators ean be expanded as a series of loeal operators. This ean be expressed as 

Oi {x+\y)Oj{x-\y)=Y, ^Ok (>') (^) ■ ( 4 ) 

k 

Here, we have used the abbreviations (x) = (vq, x) and (j) = (jo? u) for tho four-dimensional vee- 
tors. Woi^{y) are the Wilson-coefficients, whieh only depend on the relative time and distance y of 
the two original operators. The operators on the right-hand side of Eq.(4) are ordered according to 
their scaling dimensions in ascending order. This expansion works for small time differences 
(or small distances), as the Wilson coefficients behave as and be¬ 

cause the operators with larger scaling dimensions are thus suppressed by higher powers of ^|yo| 
(|y |). Eourier transforming Eq.(4), the above statement is translated into energy-momentum space, 
where the OPE is a good approximation in the large energy or momentum limit as operators with 
larger scaling dimensions are suppressed by higher powers of 1/ ^J\kQ\ (l/|fc|). 

Eor the above expansion to work in the context of a non-relativistic atomic gas, certain condi¬ 
tions have to be satisfied. Eirstly, it is important that the potential range vq of the atomic interaction 
is much smaller than all other length scales of the system, so that the detailed form of the inter¬ 
action becomes irrelevant. Eurthermore, the energy or momentum scale at which the system is 
probed needs to be much larger than the corresponding typical scales of the system. Hence, for 
the OPE to be a useful expansion, the following separation of scales must hold, which must be 
satisfied by either 1 / ^J\kQ\ or 1 /1fc|: 

rQ<^\/\/\h\, \/\k\<^\a\,na''^'^,X t- (5) 

_l/3 

Here, a is the s-wave scattering length between spin-up and -down fermions, the mean 

interparticle distance of both fermionic species, and Ar ~ 1 /\/mT the thermal de Broglie wave 
length. In other words, y/[^ or |fc| have to be large enough so that for example an expansion in 
1 / {a^J\kQ\), \/\ko\ and 1 / (Ary/|ko|) is valid, while they should be still small enough not to 

probe the actual structure of the individual atoms. 
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In practice, we will in this work take the zero-range limit ro —?■ 0, study the system at vanishing 
temperature T = 0 and will in the eourse of the derivation of the sum rules take the unitary limit 
a —oo. Furthermore, for studying the detailed momentum dependenee of the speetral-density, we 
will in the following diseussion make use of an expansion in 1/^|ko|, but not in l/|fc|. |fc| will 
instead always be kept at the order of Fermi-momentum of the studied system. 


B. The OPE of the single-particle Green’s function for general values of a 


In this paper, we will employ the OPE of the single-partiele Green’s funetion, whieh was eom- 
puted in [25]. Let us here briefly reeapitulate this result and diseuss its form rewritten as an 
expansion of the self-energy E^(ko, k). The starting point of the calculation is 


ig^{k) ^ / jye'^(r[vAt(x+|)v4(x-|)]) = 


( 6 ) 


ko-Cfc-^tW’ 

where k should be understood as (k) = {ko,k). The OPE for g^{k) can then be earned out, as 
diseussed in detail in [25]. If translational and rotational invarianee holds, all sorts of eurrents 
vanish and the OPE expression (taking into account terms up to momentum dimension 5) can be 
simplified as follows: 


^ _ G^{k)A{k)n^ ■ 


C 


C 


Anma 


{k)—A{k) - ^G\k)r.^\k,0-k,Q) 

oko I 


3 ... mj., 3^ 


Here, G{k) is the free fermion propagator, 

G{k) = 


[., ^ ■ 

(2;r)3 2m 


(7) 


1 


ko — £fe ’ 

A{k) represents the two-body seattering amplitude between spin-up and -down fermions, 

4;r 1 


A{k) = 


m 


4— mko —\/a 


( 8 ) 


(9) 


and {k^p-,k\p') stands for the regularized three-body scattering amplitude of a spin-up fermion 
with initial (final) momentum k {k') and a dimer with initial (final) momentum p {p'). “regularized” 
means that infrared divergences originally appearing in the three-body scattering amplitude have 
been subtracted (see Sections III C and III E of [25]): 


r;‘=®(k,0;k,0) = rt(k,0;k,0)-A(k) J 


( 10 ) 
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Furthermore, Po{q) is the momentum distribution function of spin-c fermions, the density of 
spin-down fermions and C the so-called contact density [28-30]. 

Comparing Eq.(7) with the definition of the self-energy of Eq.(6), one can easily find an ex¬ 
pression for which (again up to terms with momentum dimension 5) is consistent with the 

OPE of the single-particle Green’s function: 

C 

Anma dko 


Z^^^{k)=-A{k)n^ 


^ ■A{k) - ^Tl^^{k,Q-k,Qi) 

I 


r d 


0 

3 ^2 


■m+^L^Aik) 


2 r 


dq q 
(2;r)^ 2m 


pM)-^ 


( 11 ) 


-dko, ^ ^ dkf 

Assuming the considered system to be spin symmetric [pf(q) = pi{q)], the integral of the mo¬ 
mentum distribution function appearing in the above equation can be evaluated by one of the 
Tan-relations [28-30], 


E 

a=U‘ 


dq q^ 


Pa(9)-^ 


(2;r)^ 2m l 

where 8 is the energy density of the system. We hence get, 
.OPE'- ■ P- 


= 8 + 


C 


Anma 


( 12 ) 


Lp(*) =-A(*)„* - 


1 


-A{k) + 


-0 

3 ^2 

i:4^w 


8 + 


Anma/ 


(13) 


2 Vdko '' ' ^ dkj 

Among the various terms appearing in Eq.(13), the most involved piece to evaluate is the three- 
body scattering amplitude T^^^{k,0-,k,0), which will be studied next in a separate subsection. 


C. Three-body scattering amplitude 


The difficulty in obtaining 0;/:, 0) stems from the fact that this scattering amplitude by 

itself does not solve a closed integral equation and therefore can not be computed directly. We 
thus have to use T^{k,0\p,k — p) with a more general momentum dependence, which will, for 
simplicity of notation, from now on be denoted as T^^{k-,p). T^^{k-,p) satisfies the following integral 
equation (note that we for the moment work with the non-regularized version of the amplitude): 

'dqodq 


Tt{k;p) =G{-p) + i / T^{k\q)G{q)A{k-q)G{k-p-q) 


(2ny 


1 


Po + £p 

f dq 


( 14 ) 


4;r 


T^{k\eq,q) 


(2^)^ 2 3q^ — 2q k + k^ — Amiko — 2 


{p+q-kp 


+ mpo + ^ - mko 
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In going to the second and third lines, the integral over qQ is performed and thus qQ is fixed to £q. 


Next, setting po = Ep provides a closed equation, 
f dq Alt 


m 


T^{k\ep,p) = —2 


T\{k\£q,q) 


(15) 


P" J (2^)3 i V3g2 _ 2g ■ fc + fc2 - AmkQ - i (P+g-fe)^ + g!+g! _ ^k^'' 

which needs to be solved numerically. The technical details of this step are presented in Appendix 
A. Once the above equation is solved and Ep^p) has hence been obtained, one can extract the 
desired amplitude from Eq.(14) by setting p = k: 

I f dq 4;r T^{t,£q,q) 


T^{k;k) = 


ko + £k J (2^)^ 4a/ 3q2 — 2qkA-k^ — Amko — 1- 


1 


■ + 


dq 


Alt 


m 


a 


ko + Ek J (2^)^ ^a/ 3q^ — 2q k + k^ — Amko 
f dq Ak T^ik^Eq.q) A- 


(16) 


m 

7:1 


J ( 2 ^)^ 5 V-Iq-k + k^ — Amko — 2 

Finally, returning to the regularized scattering amplitude T^^^{k,0',k,0) = T^^^{k-,k) [defined in 
Eq.(lO)], we get. 


=mk-,k)-A{k) I 


m 


dq 

{2ny \q^ 


1 


■ + 


f dq 

4;r 

A7t 

1 (2;r)3 

5 a / 3q^ — 2qk + k^ — Amko — ^ 

k'^ - Amko - k_ 


m 


(17) 


dq 


Alt 


T^ik^Eq^q) + 


•> ( 27 r )3 i ^ 3 g 2 _ 2 qr. _ i q^ 

D. The OPE of the single-particle Green’s function in the unitary limit 


So far, we have studied the OPE for arbitrary values of the s-wave scattering length a between 
the two spin degrees of freedom (which should however be kept large enough for the conditions 
of a valid OPE to apply). One could in principle continue with these general expressions, derive 
sum rules for nonzero a~^ values and analyze them according to our strategy outlined in the 
introduction. 

In order to provide a clear account of the proposed method, we will however not do this here 
but concentrate on the unitary limit (a —?■ 00 ), which considerably simplifies many of the equations 
needed to derive the sum rules, but already exhibits all non-trivial technical difficulties that will 
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arise in an analogous, but more involved manner when generalizing the ealeulations to nonzero 

. 

Firstly, looking at the unitary limit of the OPE result of Eq.(13), the terms proportional to 
vanish and the faetor eontaining derivatives of A(k) ean be obtained in a simple form: 



2^/^7t £k-ko 


(18) 


As for the ealeulation of the three-body seattering amplitude T^^^{k,0',k,0), the integral equa¬ 
tion of Eq.(15) is made slightly more manageable because of a vanishing term in the first 
denominator of the integrand on the right-hand side. The regularized scattering amplitude itself, 
given in Eq.(17), also simplifies as the integral appearing in its second term [see Eq.(17)] can now 
be performed analytically: 


f dq 

47t 

4;r 

m 

J (2;r)3 

2 y/ 3q^ — 2qk + k^ — 4mko 

^ — 4mko 

9^ 


1 

7t 


\/3 ^_ 3ko —£fc / 1+ \/3Vl — 2ko/£fc \ 

2ko-efc V^(£fe-2ko)3/2 -l + V3y/l-2ko/efe/ 


(19) 


Eor a spin-symmetric system, making use of the equations of motion and Tan-relations, it is 
possible to express the expectation values of the local operators appearing in the OPE in terms of 
particle density n^, energy density £, and contact density C [see Eq.(13)]. In the unitary limit, these 
quantities only depend on one single scale, which determines the properties of the system. Here, 
we define the Eermi momentum and the Eermi energy by = n^ = k^/{6n^) and £p = kp/(2m). 
At infinite scattering length a ^ (and zero temperature T = 0), 8 and C are then given as 




% 

lOTZ^m' 



( 20 ) 


These values have by now been extracted from both theoretical quantum Monte-Carlo simulations 
and experimental measurements, which give consistent results, as shown in Table I. In the spe¬ 
cific analyses presented in this paper, we will use the values obtained from quantum Monte-Carlo 
studies (denoted as “simulation” in Table I). 

Assembling all the results and definitions of the last two subsections, we reach the following 
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TABLE I. Numerical values of the Bertsch parameter ^ and the dimensionless contact density in the 
unitary limit at zero temperature. The column “simulation” gives numbers extracted from quantum Monte- 


Carlo simulations, while the column “experiment” contains values from ultracold-atom experiments. 



simulation 

experiment 


0.372(5) [34] 

0.370(5)(8) [31,32] 

c 

3.40(1) [35] 

3.33(7) [33] 


final form for the OPE in the unitary limit, 




8 3/2 


1 


1 


+ ■ 


V£k-2ko 37t^ 
3ko — Cfe 




1 


73 


1 


ko + £k ^ 2^0 ~ £k 


^ \/^i£k-3koy/^ 

£k — ko 


log 


1 + 7371 - 2 fco/£fc \ ^ 

-l + 737l-2VeJ 


( 21 ) 


e 5/2 _ 

5^ ^ {£k—2koy/^ 


+ 0{kQ ^), 


where we, for simplieity of notation, have introdueed the funetion L{x), which is defined as: 

dq An 




mk',£q,q) + ^ 


{2ny 


( 22 ) 


2\/ ^q^— 2q-kA-k^— AmkQ 
Note that we here have made use of the fact that L{x) is dimensionless and hence can only depend 
on the ratio ^o/£fe- As mentioned earlier, L{x) can be obtained by solving Eq.(15) and substituting 
the result into the above definition. The detailed steps of this procedure are given in Appendix A. 
Here, we simply note that the imaginary part of L{x) (which is its only piece that will play a role 
in the sum rules to be derived later) is a finite, but sharply peaked function, which is non-zero only 
in the interval: 1/3 < x < 1 (see Eig. 7). 


E. Derivation of the sum rules 


We now derive the sum rules from the OPE of Eq.(21). Eor doing this, we consider ko to be a 
complex variable and study the contour integral, 

[ dko [E|(ko, k) - E5^P^(ko, k)] lC{ko) = 0. (23) 

J Cl +C 2 f 

Here, E|(ko, k) is the exact (and at this moment unknown) self-energy, E|^^(ko, k) is its approx¬ 
imate OPE expression of Eq.(21). /C(ko) is assumed to be an analytic function on the upper and 
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FIG. 2. The contours C\ and C 2 on the complex plane of Uq, used for deriving the sum rules. The wavy line 
on the real axis represents possible locations of non-analytic poles or cuts of £^(^0)^) and 


lower half of the eomplex plane of and to be real on the real axis, but is otherwise eompletely 
arbitrary. The eontours Ci and C 2 are shown in Fig. 2, in whieh the wavy line depiets possible 
non-analytie poles or outs of k) and k), whose aotual locations depend on the cho¬ 

sen value of |fc|. The above integral vanishes because the exact self-energy E|(fco, k) and its OPE 
counterpart are analytic in the upper and lower half of the complex plane. Furthermore, we know 
that the OPE is valid at large \ko\, from which follows that the integrand on the left-hand side of 
Eq.(23) vanishes (to the order we are considering) along the large half-circles in C\ and C 2 . As we 
have assumed /C(ko) to be real on the real axis, it is noted that the added contour sections along 
the real axis leave just the imaginary parts of the self-energies, while their real parts vanish. Thus, 
we can now write down the sum rules as 

/ oo poo 

<i(OlmE|(tt)-t-/0^,fc)/C(tt)) = / <i£0lmE®^^(£0-t-/0^, A;)/C(£0), (24) 

-CXI J —<x> 

where here and in the rest of the paper co is understood to be a real variable. The right-hand 
side of this equation can be calculated from Eq.(21), once the kernel IC(co) is specified. This last 
step, however needs some care, as some terms of Eq.(21) at first sight lead to divergences on the 
right-hand side of Eq.(24). This is for instance the case for the last term in Eq.(21), which has 
an imaginary part for kp = ® > £fe/2 and diverges as (co — efc/2)~^/2, when co approaches 
from above. This superficial divergence originates in our sloppiness of treating cuts in the above 
derivation and can be cured by taking into account all parts of the contours C\ and C 2 which run 
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along the cuts and their thresholds. The details of this procedure are given in Appendix B, where 
it is explicitly shown how all superficial divergences cancel and that hence the right-hand side of 
Eq.(24) is indeed finite. 

All this then leads us to the following form of the sum rules: 


<itt)/C(tt))ImE|(£0 -t- k) 


3n: ^ 


/ dCO2(0-£kK:\co) + —C£F —/C(^)-/C(-efe) 

/efc/2 3;r L ;r ' ' 


+ 


A p2 

.r 

b 


'■Efc/2 




d(o^J£k-2(o 6/C'(£o) - (Cfe-3ct))/C"(m) 


8 


15;r 


/ doojlco-£k 3K!'{(O) + {(O - £h)K!" {(O) 
3efc/2 L 


(25) 


For deriving this expression, we have, additionally to the assumptions mentioned earlier, assumed 
that /C((o) vanishes at to —)■ oo faster than 1 /a/w. If one wishes to use kernels which behave 
differently (as for instance in the so-called finite energy sum rules in QCD [37], see also Appendix 
C), one should go back to the OPE of Eq.(21) and rederive the corresponding sum rules. Our 
statement made above on the cancellations of superficial divergences however still holds for this 
case. 

Furthermore, in the limit ko = W ^ £fe> Eq.(21) takes a considerably simpler form, making it 
thus possible to derive the resultant sum rule with much less effort. Moreover, if one introduces 
certain assumptions of the functional form of the self-energy, one can even analytically extract 
some of its properties from the sum rules. How this can be done by making use of the finite energy 
sum rules, is demonstrated in Appendix C. While providing a simple and qualitatively correct 
picture, this approach however has the drawback of relying rather heavily on mean-field theory for 
fixing the form of the self-energy and therefore is inferior to the MEM analysis to be presented in 
the following sections, which does not need any other input besides the sum rules themselves. 


F. Choice of the kernel/C(ft)) 

As a next step, we have to fix the concrete form of the kernel lC{(o). As discussed in the 
previous sections, this kernel must be analytic on the complex plane of CO and real on the real axis. 
Furthermore, IC(co) should vanish faster than 1 /at (O —)■ co on the real axis. Obviously, these 
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restrictions still give room for an infinite number of choices. From the experience of QCD sum 
rule analyses, it is however known that a simple Gaussian centered at the origin works well for 
extracting the lowest poles of the spectral function. We will in this paper follow a similar strategy 
and use 

/C„((0,M) = n = 0,1 (26) 

as our kernel. M is usually referred to as the Borel mass in the QCD sum rule literature, which 
we will follow in this work, while in [23, 24] the symbol COo was used for this variable. M can in 
principle be freely chosen as long as the OPE converges. As will however be shown in Fig. 3, the 
OPE convergence worsens for decreasing values of M, which means that there exists some lower 
boundary of M, below which the OPE is not a valid expansion. 

As the imaginary part of the self-energy on the right-hand side of Eq.(25) extends to both 
positive and negative values of co and is in general not an even function, it is noted that using only 
the most simple kernel with n = 0 does not suffice to determine ImL |((0 + i0~^,k) as for this kernel 
all odd-function contributions automatically drop out of the sum rules. We hence need to introduce 
one more kernel which should be an odd function in CO, for which the n = 1 case in Eq.(26) seems 
to be the most natural choice. 

Eet us mention here that in the literature of QCD sum rules, other kernel choices have been 
proposed, such as a Gaussian with a variable center [38-40] or with complex Borel masses, which 
leads to an oscillating kernel [41, 42]. For this first study, we however prefer Eq.(26) because of 
its simple analytic form. 

Substituting the above kernels into Eq.(25) then gives the final form of the sum rules. 
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(28) 


where /v(j) and ^v(>') are the modified Bessel funetions of the first and seeond kind, respeetively. 
Furthermore, the funetions G'j(y) have been defined as follows: 
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(29) 


The ratios of the right-hand sides of Eqs.(27-28) and their respective leading order terms are shown 
in Fig. 3 as functions of the Borel mass M for three typical values of the momentum |fc|. 

The sum rules of Eqs.(27) and (28) look quite cumbersome, but their analytic structure becomes 
clearer if one takes the small momentum limit (C/j. —0). Using the kernel of Eq.(26) with general 
values of n, one can show that in this limit the EO term behaves as M^/2+« NNEO term as 

M 1/2+". The NEO term on the other hand can be shown to be proportional to M® = 1 for n = 0, 
while it vanishes for all other positive n values. The results for n = 0 and n = 1 are given by 
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FIG. 3. The right-hand sides of Eqs.(27) and (28), divided by their LO terms, as a function of the Borel 
mass M. The left and right plots show the cases of « = 0 and n = I, respectively. Starting from the top, 
each line shows the OPE for momenta |fc|/^F = 0, 0.6 and 1.2. Here, EO corresponds to the first line on the 
right-hand side of Eqs.(27) and (28), NLO to the second and third lines and NNLO to the fourth and fifth 
lines. The vertical arrows at the bottom of each plot indicate the lower and upper boundaries of the regions 
of M, which will be used in the MEM analysis of Section III. 


17 











































and 


(i(0/Ci((0,M)ImE-[.((0, fc) 

D 

_ + ^r(l/4)e''V^ - 

~ V 3 £f / 


(31) 
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Here, the term proportional to Cfc in the last term eomes from Taylor expanding the leading order 
terms of the first lines of Eqs.(27) and (28) in Cfe/M. The above equations should give the reader 
an idea on the behavior of the OPE at least for small |A;|. In the aetual analysis of the next section, 
we will however use the full result of Eqs.(27) and (28). 


III. MEM ANALYSIS FOR THE SPECTRAL DENSITY 

Next, we discuss the imaginary parts of the self-energies, which we have extracted numerically 
from the sum rules by using the maximum entropy method (MEM). This sort of approach for 
analyzing sum rules, was recently applied to QCD in a similar way [36] and has during the last 
few years been used to study hadrons in various environments [40, 43-47]. Eor the technical 
details of this analysis, we refer the reader to Appendix D and the references cited therein. 

A. The Borel window and the default model 

Before discussing our results, let us here at first briefly explain how to determine the lower 
and upper boundaries of the Borel mass M used in the analysis. Eor fixing the lower boundary 
Mmin, we demand that the highest order (NNEO) OPE term, which is proportional to ^, should be 
smaller than 10% of the leading order term. Note, that this condition generally leads to a value of 
Mmin, which depends on the momentum \k\. We will here first fix A /min at \k\ =0 and take this 
momentum dependence into account only if it leads to an increasing value of M min . This keeps the 
momentum dependence of M^m to a minimum and at the same time ensures that for any value of 
I fc I , only Borel mass ranges with a satisfactory OPE convergence are used as input for the MEM 
analysis. 

Eor fixing the upper boundary Mmax , we do not have such a clear-cut criterion and therefore 
can in principle choose it freely as long as it lies above M mm - Eor the analysis presented in this 
paper, we will set it as M^ax = +x, with x = 5£f. We have checked that our results do not 

much depend on this choice and the exact value of x hence does not play any important role in the 
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present analysis. The specifie values of A/min and Mmax for some typieal momentum values are 
indieated in Fig. 3 as vertieal arrows at the bottom of each plot. 

As for the default model m{(o), which is an input of the MEM algorithm (see Appendix D for 
details), we will use 


/ N 4\/2 3/2 


1 


(32) 


((02+j)1/4’ 

with y = gp. As can be understood from Eq.(21), the above default model approaches the correct 
asymptotic limit of ImL|(( 0 ,fc) ~ —{4\/2£^^^)/{3Ky/a), as ft) ^ gfc and is therefore a suitable 
choice for the present analysis. Eor avoiding singularities at ft) = 0, we have introduced the pa¬ 
rameter y for smoothing out the function around the origin. We have tested different choices for y 
and found that this affects our analysis results only very weakly. 


B. The single-particle spectral density 


After these preparations, we can now finally proceed to our analysis results. Eirst, we show the 
imaginary part of the self-energy, for three representative momenta in the left column of Eig. 4. Eor 
illustration, we show in these plots also the used default model of Eq.(32). It is seen that for zero 
momentum, the spectral function is composed of one single peak around ft) = 0 and a continuum 
behaving as ~ 1/ y/a in the positive energy region. As the momentum increases, the initial peak 
separates into two distinct peaks which start to move into opposite directions. The continuum also 
recedes into the positive ft) region with increasing momentum, leaving a growing region around 
the origin without any strength at all. 

With the extracted ImE|(ft), k), we next compute the real part of the self-energy by using the 
Kramers-Kronig relation 


/ 1 f°° TmEt(ft)',fc) 

ReEt ft), k) = -P / Jft)'-, 

' 7t J-oo ft) - ft)' 


(33) 


and executing the principal value integral numerically. The result of this evaluation is given in the 
middle column of Eig. 4, where we also show the curve ft) — e^, which appears in the denominator 
of the right-hand side of Eq.(3). It is clear from this equation that if the imaginary part of the self¬ 
energy happens to be small, the single-particle spectral density will have a narrow peak wherever 
ReE/(ft), k) coincides with ft) — e^. 

As a last step, we simply plug the real and imaginary parts of the self-energy into 


A|(ft), k) 
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1 

ft) 4- )0+ - Cfe - E|(ft) -1- i0+,k) ’ 


(34) 
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FIG. 4. Left column: Results of the MEM analysis of Eqs.(27) and (28) are shown as red lines, while the 
used default model [see Eq.(32)] is indicated in blue. Middle column: The real parts of the self-energies 
obtained from Eq.(33) and are plotted as red lines, and the function ft) — is given in blue. 

Right column: The spectral density A-^{( 0 ,k), as computed from the results of the two columns on the left 
and Eq.(3). As in Eig. 3, each row from top to bottom corresponds to momenta |fc|/^F = 0.0, 0.6 and 1.2, 
respectively. 


to obtain the single-particle spectral density A^{(o,k). The resulting functions are given in the 
right column of Fig. 4. It can be seen there, that for small momenta \k\, the spectral density 
is dominated by the narrow hole-branch in the negative energy region, while the particle-branch 
consists of only a relatively broad bump. This changes at around |fc| ~ 0.5fcp, where the main 
strength of the spectral density switches over to the particle branch, which, as the momentum is 
further increased, starts to move into the positive energy direction. On the other hand, the hole- 
branch bends back into the negative energy region, while gradually losing its strength. To give the 
reader a better visual grasp of the spectral density as a whole and especially on the behavior of 
the particle and hole branches, we show A^(( 0 , k) in a density plot as a function of both energy ft) 
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FIG. 5. Density plot of the spectral density A-^{( 0 ,k) shown as a function of energy ft) and momentum \k\. 
The green dashed lines indicate the results of a fit of the particle and hole peak-maxima to Eq.(35). 

and momentum |fc| in Fig. 5. To improve the visibility of this plot without ehanging its essential 
features, we have artifieially inereased the imaginary part of E|(ft), k) in Eq.(34) by an amount of 


0.2 £p. 


In this figure, the typieal BCS-like dispersion of the partiele and hole branehes elearly manifest 
themselves. Qualitatively, this result agrees with the speetral densities extraeted from both quan¬ 
tum Monte-Carlo ealeulations [7] and a Luttinger-Ward approaeh [6]. In order to make a quantita¬ 
tive eomparison with other methods, we fit the peak maxima to a dispersion relation parametrized 
as 



(35) 


whieh we have adopted from [6]. The resultant eurves are shown in Fig. 5 as green dashed lines. 


while the eorresponding values of /i. A, and are given in Table II. It is seen in Figure 5 


that the fit is able to reproduee our dispersion relation fairly well, with the exeeption of the low 
momentum region of the partiele braneh, whose eurvature ean not be eaptured fully by the simple 
formula of Eq.(35). Note that this leads to a slight overestimation of the gap A. If we simply read it 
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TABLE II. Fit results of the particle and hole branches shown in Fig. 5 to a dispersion relation parametrized 
as in Eq.(35). 





Particle 


Hole 


fi/ep 

A/ep 

m+/m 

U+/ev 

m jm 

U-lev 

this work 

-0.18 

0.57 

1.02 

-0.37 

1.09 

-0.12 

[6] 

0.36 

0.46 

1.00 

-0.50 

1.19 

-0.35 


of from the point at which the particle and hole branches are closest, we get a value of A/Cp = 0.54 
instead of the one given in Table II. 

Comparing the values of this work with those of [6], it is seen that the two approaches give 
comparable results for the gap parameter A, effective masses and Hartree shifts for both 
the particle and hole branches. On the other hand, the chemical potential /i deviates significantly 
from [6], even giving a different sign. The reason for this discrepancy apparently originates in the 
low sensitivity of the sum rules to the absolute position of the CO axis. This can be understood by 
inspecting the OPE of Eq.(21). After setting ko = CO and making a change of variables CO ^ Co' 
as CO = to' + cop, with cop of the order of Cp and expanding the resulting expression in coq/co', one 
notes that only the NNEO term of the OPE will be modified, which must be kept small due to the 
convergence condition of the OPE. Therefore, we can expect that such a change of variables will 
introduce no qualitative modification of the OPE, while the spectral density experiences a parallel 
shift of cop. 

It is in principle possible to choose cop such that the fitted value of /i approaches the correct 
value of around 0.36 Cp. Due to the convergence criterion of the OPE, such a choice however 
leads to a significantly larger value of Mujin and therefore to a rather poor resolution of the MEM 
extraction of ImE|(co, k). We have thus not explored this possibility any further and simply note 
that at the present stage, the absolute positions of the structures appearing in the spectral density 
should not be taken too seriously. 

As a final point, we study the density of states of the single argument co, P|(co), which is 
obtained by integrating the spectral density over the momentum | fc |: 

f (Pk 

Ptio}) = J j^A^{co,k). (36) 

This function is shown in Eig. 6, from which one can immediately read off the approximate gap 
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FIG. 6. The density of states, Pf (ft)), obtained by integrating the spectral density A-f{(0, k) over the momen¬ 
tum k as shown in Eq.(36). 


value, which can be regarded as half the width of the region where P|(ft)) loses almost all of its 
strength. To draw Fig. 6, we have added an constant amount of 0.002 Cp to the imaginary part 
of E|(ft), k), which reduces artificial effects caused by evaluating the integral numerically from a 
discrete number of data points, but does not change the gap structure of this plot. 


IV. SUMMARY AND CONCLUSION 

The work presented in this paper was carried out with two essential goals in mind. As the 
introduced techniques are new and have not been applied to cold atom systems so far, we first 
needed to test to what extent the sum rules and MEM are able to extract the single-particle spectral 
density from the result of the OPE. This is by no means a trivial test, because the OPE considered 
in momentum space does not converge for momenta below the Fermi momentum [25], as we have 
already discussed in the introduction. It was therefore at the beginning not clear to what degree the 
sum rules can extend the applicability of the OPE to lower momenta or energies. As it however 
turns out, even at zero momentum |fc| and small ft), the sum rules of Eqs.(27) and (28) lead a fairly 
reasonable behavior for the spectral density, which suggests that our approach is indeed useful for 
extracting the spectral density at any momentum and energy. 
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Once the proposed method is proven to work well, our second goal was to provide an indepen¬ 
dent framework for evaluating the superfluid pairing gap A of the unitary Fermi gas. Our obtained 
value is given in Table II and can be inferred from Fig. 5. We wish to emphasize here that even 
though we have only taken into aceount the first few terms of the OPE, in whieh the Bertseh 
parameter and the contact density are the only input values, our numerical result shows reason¬ 
able agreement with other theoretieal approaehes [6, 8]. Speeifieally, we obtain A/gp = 0.54, 
when extraeting the gap from point of smallest distance between the partiele and hole branehes 
and A/gp = 0.57 from an overall fit of our dispersion relation to Eq.(35), while [6] and [8] get 
A/gp = 0.46 and A/gp = 0.50(3), respeetively. Eor eonfirming these results in the future, it will be 
neeessary to eonsider still higher order terms in the OPE, evaluate the size of their eontributions 
and examine their impaet on the speetral density. 

Using the method proposed in this work, we have so far only studied the fermionie single- 
partiele ehannel at zero temperature. As long as the eonditions for its applieability (that is tq 
1/\/|^o| , Xt) are satisfied, the OPE technique is fairly general and can in principle 

be applied to any kind of bosonie or fermionie systems with one or more constituents. One ean 
therefore envisage various future applications of this approach. Eor instance, in [23] the OPE for 
the retarded eorrelator of the density operator has already been worked out, and one in principle 
just needs to apply MEM or some other sort of fitting method to extract information on the dynamie 
strueture faetor from the OPE expression. Another interesting direetion of research could be the 
generalization of this approach to finite temperature. Eor being able to do this, one however 
needs information on the finite temperature behavior of the operator expeetation values whieh 
appear in the OPE of the ehannel of interest. Eor the system considered in this paper, this would 
eorrespond to the finite temperature values of the Bertseh parameter and the eontaet density, whieh 
are ealeulable using quantum Monte-Carlo simulations [48]. 
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Appendix A: Numerical solution of T^^^{k,0-,k,0) in the unitary limit 


As discussed in Section IIB, we need to solve the integral equation given in Eq.(15) numeri¬ 
cally, which is then substituted in Eq.(17) to obtain the desired scattering amplitude 0; k, 0). 

The technical details necessary for this task will be outlined in this Appendix, in which we gen¬ 
eralize the discussion given in Appendix B of [25], where ko was set to £k, while we here have to 
keep it as an independent variable and will specifically consider ko = CO + iO^, with CO being real. 

Eirstly, it is noticed that the dimensionless function 


s^{k;£p,p) = — T-^{k;£p,p+^k) 


(Al) 


satisfies a simpler integral equation, which is given as 


st{t,£p,p) 


(p + ifc)2 

f dq 


16;r 


^{k\£q,q) 


(2^)^ J^q^ + lk'^-Amko V + 2q2 + 2p ■ q -h ifc2 - 2mko 


(A2) 


= -I{k;p) 


dq 


■J{k;p,q)s^{k;£q,q). 


J (2;r)3' 

The important point here is that the Kernel J{k-,p^ q) now depends only on the angle between q 
and p, which will permit a partial wave expansion of the above integral equation. 

Next, we expand s^{k-, £p,p) into its partial waves, which depend on the angle Q between k and 
p as 

oo 

s^{k-,£p,p) = Y,s^H^MPi (cos e ), (A3) 


1=0 


where (x) are the Eegendre polynomials. We have made use of the fact that is a dimensionless 
function, which can hence only depend on the ratios ko/£fe and £p/£k- It can be shown that each 
term {ko/£k,£p/ Cfc) in the sum of Eq.(A3) satisfies a closed integral equation. 


„(0 f h _ f 'jd) ^ £e £ 2 . 

t \£k^£kj \£kj Jq \k\'^ yefc’efc’efc 


„(0 f kp £q 

t \£k^ £k)' 


(A4) 


Here, the function is defined as 

J ^d cose Pi {cos e)X{k-,p), 

which can be rewritten with the help of the Gaussian hypergeometric function 2 X 1 («, t>;c;y): 


(A5) 



l\ 9 
{2l-l)\\l+9x 


bi/T 

l+9x 



Z + l / + 2 3 36x 

^’^’^'^2’(l+9x)2 


(A6) 


25 














Furthermore, we have defined as shown below: 


'j(l) f hi ^ A 
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Using again the Gaussian hypergeometric function this gives, 
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(A 8 ) 


As a next step, we need to solve Eq.(A4) numerically for general values of 1. In practice, we 
however will not deal with this equation directly, but first define 



_ „(0 ( h) 
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(AlO) 


and then solve this equation for 5s^^\kQ/ek,,£p/£k)- This is done in order to avoid (or at least to 
weaken) the singularities that appear in ^^(ko/SfcjSp/Sfe) for certain values of ko/£fe and Cp/Cfc 
and use instead the better behaved 5s^^\kQ/£k,£p/£k)- Once this is done, the result is substituted 
into Eq.(22), which, by making use of the above definitions, can be rephrased as 


(-) 

\£kJ 


= £fe 


dq 


Stt 


rt(ke„q) + p 


(27r)3 3 / 3q^ — 2q ■ k + — 4mko 
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d 
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3£q I i _ ^ 
2 efc ^ 3 Efc 


\ J T Y 


(All) 


After obtaining 5s^^\kQ/£kT£p/£k) for each value of I individually, the corresponding contri¬ 
butions are added in Eq.(Al 1), which then gives the final form of L(ko/£fe)- It suffices to evaluate 
the function L(ko/£fe) for one specific value of e^, as its form for general can be obtained by a 
simple rescaling of its argument. 

In Eig. 7, we show the final results for Im \L{x)\ for various maximum values of I in the sum of 
Eq.(Al 1). (We show only the imaginary part because this is the only piece that will be needed for 
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FIG. 7. Im [L(;c)] , obtained by solving Eq.(AlO) and adding up the results in Eq.(Al 1) to various maximum 
values of /. 

constructing the sum rules.) It is seen in these plots that Im [L(x)] is essentially determined by the 
first 5 terms in the sum over I and that the expansion eonverges quiekly for values beyond / ~ 10. 
In this work, we will take into account the terms up to / = 20. 

Furthermore, it ean be seen in Fig. 7 that Im[L(x)] takes non-zero values only in the interval 
^ < X < 1, where it peaks sharply at around x ~ 0.5. 


Appendix B: Derivation of the sum rules for a generic kernel 


To derive the general form of the sum rule, given in Eq.(25), we need to compute the right-hand 
side of Eq.(24), or, to be more precise, need to evaluate the contour integrals of k) 

along the seetions of the eontours C\ and C 2 , which run above and below the real axis. The OPE 
expression for the self-energy is given in Eq.(21) of Section IID and is reproduced here once more: 




3/2 
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VCfe —2ko 
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ko + fife ^ 2ko — fife 


^ \/^(£fe -2ko)3/2 

8 e 5/2 fife - ^0 




(£fc -2ko)^/^ 


+ 0{kQ ). 


(Bl) 


The kernel /C(ko) is assumed to be analytic in the whole eomplex kq plane and to vanish faster that 
1 / y/kQ at ko —> on the positive real axis. 
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As some of the derivations are somewhat involved, we will consider each term of the OPE indi¬ 
vidually. As in the main text, we here consider ko to be a complex variable, while CO is understood 
to be purely real. 


1. Leading order (LO) 


Using 


we immediately get. 


Im 


^/£^,-2co-iO+ 


= e{2(o-£k) 


V2(0 - £k ’ 


dco)C{co)lmL\^{co + k) 
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■Jek/2 


„„ -=/C((0) 

e^/1 ^2(0 - £k 

D _ 

dco^/2co-£kJC'{co). 


(B2) 


(B3) 


Note that for the above integrals to converge, the assumption of /C((o) to approach 0 quicker than 
1 /^/a at CO —oo is needed here. 


2. Next-to-leading order (NLO) 

Being proportional to the contact density parameter the NLO expression consists of two 
pole terms, one log-term and one term containing the function L{ko /fife). The pole terms are easily 
treated using 


Im 


1 


CO-x + iO+ 


= -7td{co-x), 


(B4) 


which gives 


<itt)/C((0)ImE^^^’^°'®((0 -h zO^', k) 




(B5) 




3;r L 2k 

Next, we consider the log-term, which needs a somewhat more careful treatment of the contour 
integral, because simply taking its imaginary part leads to a divergence at £0 = efc/2. Before doing 
this, we note that 

1 -I- a/3^/I —'lcol£i. 

(B6) 
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— 1 -(- \/3 1 — 2co j £k 
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FIG. 8. The contour integral on the complex plane of the variable needed to calculate the NLO log-term 
contribution to the sum rule. 


which is the argument of the log in Eq.(21), is positive and real for (0 < Ef./3 and negative and 
real for e^/S < co < Ekj2, where the log therefore has a eut. On the other hand, for co> Ej. 12, the 
above expression ean be rewritten as follows: 


f + ^/3^/Y^^2coJ£k^^l^ 1 -z\/3A/2(o/£fc- 1 

-1 -f-i^/3^/2(oJ£k^ 

= 3m-£fc 


where 0 is given as 


0 = tan ^ 


/ ^3£^y/2C0 - Ek \ 

V 3co-2Ek /' 


(B 8 ) 


Therefore, the log of Eq.(B 6 ) is purely imaginary for ft) > £fe/2. In this region, the root in front of 
the log in Eq.(21) is also purely imaginary, whieh means that the term as a whole is real and that 
there is no out for co 

Henoe, it is understood that we just have to evaluate the oontour shown in Eig. 8 . The oorre- 
sponding analytioal formula is 


JC\ —C 4 


2 / JCi~C4. 


dkoKiko) 


3ko - £fc 
{Ek-2ko)'^/^ 


log 


l + y3v/l-2fco/£fc \ 

-l + x/3v/l-2Wefey’ 


(B9) 


for whieh we below ealeulate the parts C\ - C 4 separately. 
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Firstly, it is seen that the integrand is not singular at = £fe/3. Thus, the eontour Ci eireling 
around this point vanishes as its radius approaehes zero: 


^ jf dko}C{ko)T^^^^^°^{ko,k)=0. 


(BIO) 


Next, the eontour segments C 2 and C 4 are eonsidered. They have a finite value due to the eut of 
the log, whieh ean be evaluated as follows: 


yj dko]C{ko)L^^^’'°^iko,k) 
JCo-\-Ca 


) C2+C4 

4 

^ tp 
I s 


'■efc/ 2 -e 


\/£fe J efc/3 

1 


3;r2 

1 




dcoIC{co 
1 


3(0 —y^ 


2v/2 


(y^ — 2 (o)3/^ 

Viicif) 


(BID 


'■Efc/2 


4e;;,/3 


(fco^/efe-2(o 6/C'(co) - (Cfc-3(o)/C"((o) 


Here, e stands for the radius of the eirele around = £fe/2 of the eontour C 3 . The last eontribution 
eomes from C 3 , whieh, after a ehange of variables from ko to 0 (fco = £fe/2 + fie'®), is divided into 
two parts: 

^ jf dcoJC{ko)I^^'^'^°^{ko,k) 

= 1 rj0e'®/C(efe/2 + ee'®)Ej^^°’^°®(efe/2 + ee'®,A:) (B12) 

2 Jk 

+ ^ [ ""dee‘^lC{£k/2 + £e‘^)L^^^’^°^{£k/2 + £e‘^,k). 

2 Jo I 


The first part is evaluated as 


^ / dee^^]C{£k/2 + £e^^)L^^^^^°^{£k/2 + £e^^,k) 
2 Jk ' 


:C 


fO 


37r^ \/£fe J7t 




xlog 


2 /^+ \/3v2ee^ 


,-V^+V3V2ee'(^-^) 


= ^C4V<iK(f )(!-/■ ‘ 


1 


y/£ \/2%^ 




CeF/C(f), 


(B13) 
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while the second one gives 


2 Jo 


' dQe^^Kie^jl + ee'®)E^^°’*°®(efe/2 + ee'^ fc) 




2 ^ £i 


Stt^ \/^fc J7t 


^ 2 eeKe+^)]^/^ 


xlog 


1 /^ + \/3\/ 2£e'(®+^) 
-\/^ + \/3\/ 2ee'(®+^) 




(B14) 


Adding the two results from above, we finally get 




2 / Jci 


\/£ y/STt'^ 


Thus, assembling all the contributions, we can obtain the result for the whole contour of Fig. 8: 


2/ dko}C{ko)L^^'^''°^{ko,k) 

J Cl —C4 


\/3n^ 


CeF[ 2 /C(f)-/C(f) 


(B16) 


+ 


rC 


'■Efc /2 


StT^ \/£fc J £fc/3 


d(0\/£k-2(0 6/C'((o) - (Cfc-3 £o)/C"(£o) 


Note that all divergences have vanished in this final expression. 

The last term that has to be considered contains the function L{ko/£/.). As its imaginary part 
has no divergences, it is straightforward to evaluate the corresponding contribution, we just need 
to take the imaginary part of L(ko/£fc) (shown in Fig. 7), multiply the kernel and numerically 
perform the integral over ft): 


[ Jft)/C(ft))Imi: 5 ^^°’'^^^"/''^^(ft) + / 0 +,fc) 

J — 00 


^ r 

37t^^ £k J-o 

J_ri 

3n^ ^ £k . 


<ift)/C(ft))Im 


r£k 


^kfi 


dcolC{co)lm\L{^) 


(B17) 


In the last line we made use of the fact that Im[L(j£:)] only has non-zero values in the region of 

^ < x < 1. 
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FIG. 9. The contour integral on the complex plane of the variable needed to calculate the NNLO contri¬ 
bution to the sum rule. 


Together with the pole- and log-terms, we henee ean eolleet the NLO results as follows: 




+ 


4 ej 

b 


'■efc/2 


\/£fe ■4 efc/3 


do)^/£k-2(0 6/C'((o) - (Cfe-3(o)/C"((o) 


(B18) 




Let us here briefly draw the attention of the reader to the faet that the term eontaining /C(efe/2), 
whieh appears in both the pole term result of Eq.(B5) and the expression of Eq.(B16), happens to 
exaetly eaneel and does therefore not show up in Eq.(B18). 


3. Next-to-next-to-leading order (NNLO) 


As in the last seetion, we here again have to eompute a eontour integral in the eomplex plane 
of ko- This eontour is shown in Eig. 9. We henee have to ealeulate the following integral: 


yj dko)C{ko)L^^^^iko,k) 

JC\ —C 3 


A dcoJC{co) 


1 


£f, — CO 


(B19) 


,,.r. , '"'(£t/2-(0)5/2 ■ 


5;r ^ 2 z/ci-C3 

Eirst, the eontours of Ci and C 3 are eonsidered. Added together, they reeeive only eontributions 
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from the cut of the root-function. This then leads to 


1 


2 / /C 1 +C 3 

72 




(^Cp ^ / d(olC{(o) 


Cfc — (0 


'efc/ 2 +e 


£0 






15;r 

8 

15;r 


1 


_+ ^££5/2 


2/C(f)-£fe/C'(f) 


7^ 


(^£p'^^ / d(0\/2(0-ek 31C"{(o) + {(o - ek)K!"{(o) 

Jek/2 L 


(B20) 


where, in similarity to the last subsection, £ stands for the radius of the circle around Ct) = £fe/2 of 
the contour C 2 . 

Next, the contribution of C 2 is calculated, leading to 




2 / /c 2 
10 ;r 


15;r 




^£y"£t/c(f) 


5 / 2 , 


-2tz 


d9e'^}C{f + ee'^) 


ie\ — 


iO 


(£g!( 6 +?r)p /2 


(B21) 


1 72 g 5/2 


2 /C(f)-£fc/C'(f) 




Therefore, the final form of the NNLO contribution to the sum rule is found to be 


dko]Ciko)L^^^^{ko,k) 

J Cl —C3 

O , PCK> _ p 

dco^yico-Ek 3/C"((o)+ ((0-£fe)/C'"((o) 


15;r 




Jek/2 


(B22) 


where, as before, only the finite term remains, while all other divergent contributions at £ —;■ 0 
cancel. 
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4. Summary 


Collecting all the terms of the last few subsections, we get the following form of the sum rule 
for a kernel /C((o), which at tt) —)■ co has to approach zero faster than 1 / ^/a: 


<i£ 0 /C(ft))ImE ^((0 + iO^,k) 


_ 8 3/2 


^oo A p \f\ 

/ da^lco - £kJC\(o) + —Cep —/C(^)-/C(-efe) 

/efc /2 jTT I TT ^ ^ ' 


4 

+ —C p 


'■efe /2 


yf £k J efc/3 
4 „ Cc 


do)y/£k-2(0 6/C'((o) - (Cfe-3(o)/C"((o) 






15;r 


£fc /2 


dC0y/2(0-£k 3JC"{(0) + {cO-£k)lC”\co) 


(B23) 


This results corresponds to Eq.(25) of the main text. 


Appendix C: Finite energy sum rules for the unitary Fermi gas 

In this Appendix, we will demonstrate how to apply the finite energy (FE) sum rule approach 
[37] to the kQ = co^ £k limit of Eq.(21). This limit will considerably simplify the analysis of the 
sum rules and, after introducing certain assumptions on the functional form of ImE|((0, fc), will 
even allow us to study them analytically. 


1. Large frequency limit 


To take the ko = ® 3> Cfe limit in Eq.(21) and expanding the result in Cfe/ft) is mostly straight¬ 
forward, the only exception being the L{ko/£k) term, which is related to the three-body scattering 
amplitude. For evaluating this term, we need to solve Eq.(15) at fc = 0 (and = 0). This integral 
equation can be rewritten as 


mko,0;£p,p) 

_ m__ f __ r|(ko,0;£q,qr) 

J ( 2 ;r )3 ^ ^3q2 _ 4,^^^ + q^+p-q-mko 

_ ^ ^^^f P^ + q^ + \p\\q\-mko '\^ 

nJo ^ \p\ y/3q^-4mko ^^\P^ + (l^-\p\\q\-mko) 


r|(ko, 0 ;eq,q). 


(Cl) 
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It can be understood from the last line above that r|(fco,0;£p,p) can only depend on |j?|. We hence 
define the dimensionless function: 

rt(fco,0;£p,p) = ^tt(IPl) (C2) 

Kq 


«-0 

and rescale the momentum in units of \/mko. The integral equation thus becomes 

/I IN 12 \q\ 1 f + \p\\q\ - 

^ ^ IpI ^3q2_4 

which numerically determines (| q |). 

The term containing the function L{x) of Eq.(22) can then be given by 

dq An ^ 


1 , / N fc = 0 


ihl) 


J (27r)3 

ySq^ —4mfco 

14 r^i I 1 

= i-/ d\q\ , 

ko 7t Jo ^ 3 q 2 - 4 


ffl 

T^{kQ,0\eq,q) + ^ 


't(l9l) + ^ 


By using the numerically obtained t|(|q|), we find 

fc=o^ 0.396797 
ko 




Together with the other terms, we hence reach the desired limit: 


-OPE 


{Jcq = ft), k) 


8 j,/2\ 1 Efc 1 

P 


-\/-2£0 2 (n/-2 (0) 

4 ^ 44 ( 1 -^- 0 . 396797 ) 

1 




4^2 3 / 2 ^ 0.20750 

^ n /—(0 ^"^2 P ^.-. 


5n P 

2. Ansatz of self-energy spectral function 


3;r2 " " £0 

_ ^£5/2 _ 5 fife \^_ ^2 

5 ;r ^ V 3 ep/ 


£0 


- 2 ' 


(C3) 


(C4) 


(C5) 


(C 6 ) 


For streamlining the notation, we first rewrite the result of Eq.(C 6 ) as follows, 

1 




^ 1 ^ 

~C 2 -h—C3- - , 

£0 2 (v^)3 


(C7) 
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where we have defined 


Cl 


8\/2 3/2 


C2 


0.207498 2 

3;r2 


C3 


2 V 2 

15;r2 


3<^ -5 


£fe \ ^ 5/2 
£f 


(C 8 ) 


As Eq.(C7) is valid at large CO, we can immediately read off the asymptotic behavior of ImE |((0 + 
iO^, k) in this region as 


Here we used 


ImE|(m + iO+,fc)~ 2^3(ym)3' 


(C9) 


Im 


1 


_V-(0-iO+ 
1 


= 0(m)- 


1 


Im 


Im 


(0 + i 0 + 
1 


.(\/-(0 -/ 0 +)' 


= -n 5 {(o), 
= -e(co) 


(CIO) 


1 




One can take the simplest ansatz for ImE |(£0 + i0^,k) satisfying the above behavior as 


ImL|(ft),fc) = 0((O —5thr) 


2 y/a 2 (a/®)3J ’ 


(Cll) 


with some parameter i'thr- (We here assume i'thr > 0.) 

However, it turns out that the finite energy sum rules for this ansatz does not have a physical 
solution for Cfc/ep < 3(^ /5. Also we already know from the BCS theory in the mean-field approxi¬ 
mation (MFA) that ImE|(£(), k) has a peak at negative CO, which is absent in Eq.(Cl 1). We are thus 
tempted to take the modified ansatz, given by a naive summation of the continuum (Cll) and the 
peak in the MEA, 


ImE|(( 0 ,fc) = -nC45{co + £k - 2^) + 9{co -sthi) 


1 ^ 


(C12) 


where C 4 = is the result in the MEA, which is also expressed using the contact density C as 


r = ^ = 

^ m2 3;r2 


(C13) 


in our units. 
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3. Derivation of finite energy sum rule 


From the imaginary part of the self-energy ImE'[-(( 0 , k), we can obtain its real part through the 
Kramers-Kronig relation 


^ ^ , ,, , ,ImE+(£(j -fzO+,fc) 

ReEf(( 0 ,fc) =-P / dco - - -. 

' 7t J-oo CO-CO' 


(C14) 


Using the ansatz of Eq.(C12), the integral in the right-hand side above reduces to 

C 4 „ r dco' 


ReE|(( 0 ,fc) = 


+ P 


+ P 


Uth. a)-co' 
r dt 


Cl 1 C3 1 

■ 2 2 


03 




(O + efe-2^ 

C 4 

(O + efe-2^ ' ^ 

(O + efe-2^ y^co V 


dt 


(C15) 


where, in the second line, we set t = \/W. The integral can be performed for (O > 0 and CO < 0, 
respectively, as 


C 4 


+ 


C3 1 


ReE|(( 0 , k) = 


_ J £0 “h £fc 2<^ \/3'thr ® 


Cl 




1 


C4 


C3 1 

® T £fc 2(^ y/Sfhi CO 


+ 


+ Cl - — 


CO ) 2^/® 

C3 


log 


\/® \/ ‘^thr 


\/w + \/‘^thr 


1 


CO J \2 


tC — 1 / ‘^thr 

-tan ' 


-CO 


>> 0 ) 
< 0). 
(C16) 


For sufficiently large CO S> 5thr> using 

\ —X 


log 


1 X 


0(1 

= - 2 a:| 1 +- + - + 


(C17) 


with a: = y/Vthr/ CO <C 1, the right-hand side of Eq.(C16) can be expanded as 


(0 V CO J yTtfeCt) V ® 


y/'^thr ^ (\/3'thr) (y/^'thr) 


.(C18) 


CO 3co^ 5co^ 

By comparing the coefficient of l/co with that in Eq.(C7), we arrive at a constraint for > 0: 


Cl (1/i'thr)^ + C-y/Vthr + C3 — 0, 


(C19) 


where 


^ ^ ^ 3.792502^ ^ 

C = C2+C4- - ^^^2 -^ 


(C20) 
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4. Solution of finite energy sum rule 


For Eq.(C19) to have a real solution for ^/s^, the following condition is necessary: 

D = c2-4CiC3 >0. (C21) 

For ^ = 0.372 and ^ = 3.40, one can check that this condition is satisfied for any e^. Then 
Eq.(C19) can be solved as 


-C±^JC^-AClC^ 

-C- Vc2- 4 C 1 C 3 

2 Ci 


(^ < C 3 < 0 ) 
(C3 > 0) 


(C22) 


where the signs are chosen such that is positive. For the smoothness of the solution of 5thr 
around C 3 = 0, we assume to take 

-C- Vc2-4CiC3 


= 


2C, 


(C23) 


for any e^. 

In summary, we find 


L-^{(O + i0^,k) 
C 4 




1 


(0 + £ k - 2 ^ y /^(0 


log 


y/ ‘^thr 


(0 + £k-2^ 

for £0 > 0 and to < 0 , respectively. 


CO J 2^/a Va + y/s^ 

1 , -1 /~^hr 

— — tan ' 


z;r 0 (tt) — i’thr) 


ft) J \2 


-CO 


(C24) 


5. Single-particle spectral function 


Now we compute the single-particle spectral function defined by 


1 


Af{co,k) =-Im 


1 


n [£ 0 -t-zO+-efc-E|((o-i-zO+,fc) 
From the expression of E|(ft) -f z0+, A;) in Eq.(C24), ^|(ft), k) reads 


{CO - COn) 


wT'i'^ft),/^^ — ^ 


ImE|(ft) + z0+,fc) 


^ [(O-Cfe-ReE|((O + z0+,fc)] -h [lmE|((O + z0+,fc 


(C25) 

^ft) < ^thr) 

(ft) > i'thr)- 

(C26) 
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CO [Ep] CO [Ep] 

FIG. 10 . Left column: Real parts of the self-energy k) [see Eq.(C 24 )]. The intersections with the red 
line (ReZ = (o — £fc) represent the peak positions of single-particle spectral densities [see Eq.(C 27 )]. Right 
column: Single-particle spectral densities A-[{(0,k) [see Eq.(C 26 )]. Eor better visibility, the delta functions 
d(co — (On) are approximated by with t = 2000 . Each row from top to bottom corresponds 


to \k\/kp = 0.0, 0.6 and 1.2. 


Here the pole(s) co = (0„(y) {n= 1,2, • • ■) are the solution(s) of 


(O-efe-E|((O + z0+,fc) =0, 


(C27) 
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FIG. 11. Peak positions of single-particle spectral densities in the (| fc |, ft)) plane (red lines). For comparison, 
mean-field dispersion relations, co = — ^)^ + C 4 , are shown by the blue lines. The thin black line 

shows the position of m = . 


and the residue(s) F„{y) are given by 


^-1 


= 1 - 


d'L^{(0 + iO^, k) 


dco 


(C28) 


co=co„ 


Figure 10 shows the plots of ReE^((0 -1- i0^,k) and k), respeetively, for ^ = 0.372 and 
^ = 3.40. The peak positions of the single-particle spectral densities as functions of |fc| are shown 
in Fig. 11, where, for comparison, we also show the mean-field dispersion relations, (0 = ^ ± 
a/(£ fe — i^)^ + C 4 . By investigating the point, at which the particle- and hole-branches approach 
each other most closely, we obtain a pairing gap value of 0.65 £f, which is not much different from 
the mean-field result, a/Q ~ 0 .68 Cf- 


Appendix D: The maximum entropy method 

Let us here briefly recapitulate the basic steps of MEM and especially explain the differences of 
our analysis to the application of MEM to statistical Monte-Carlo data. Eor more details, consult 
for instance [36, 44, 49, 50]. 

The problem to be solved with the help of MEM is given in Eqs.(27) and (28). As, however, 
the OPE on the right-hand side of these equations is only known with limited accuracy and is 
moreover only valid in a finite range of the Borel mass M, the problem of obtaining ImE|((0, k) 
from the OPE is ill-posed and cannot be solved analytically. 
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MEM now uses Bayes’ theorem, by whieh additional information about ImL|(£(),fc) sueh as 
positivity and its asymptotie behavior at large energies ean be ineorporated into the analysis and 
by whieh one then ean extraet the most probable from of ImE^((0,fc). Bayes’ theorem ean be 
expressed as 

(Di) 

where H denotes prior knowledge of ImE'[-((0, k) and /’[ImE|D//] represents the eonditional prob¬ 
ability of ImE'[-(( 0 , k) for given k) and H. Maximizing the above funetional with respeet 

to ImE|(£ 0 ,fc) will provide the most probable speetral funetion. P[D\lmLH] is ealled the “likeli¬ 
hood funetion” and is obtained as 


P[D\lmLH] = 

LilmE] = --^-- / dM ^ ^ , 

J |(A/, fc) 

with n = 0 or 1. Here, fc) is given on the right-hand sides of Eqs.(27) or (28), while 

k) is defined as 

/ oo 

J(0/C„((0,M)ImE|((0,fc), (D3) 

-oo 

and henee implieitly depends on ImE'[-((0, k). The error funetion |(M, k) stands for the uneer- 
tainty of fc) at Borel mass M and momentum |fc|, whieh we determine from the uneer- 

tainties of the parameters ^ and ^ (e.g. the Bertseh parameter and the eontaet density) appearing 
in the OPE. 

P[ImE|//] on the other hand is ealled the “prior probability” and ean be written down as follows: 


P[lmL\H] = 

poo p 

S[ImE] = / dco ImE^((0,fc)—m((o)—ImE|(ft),fc) log 


ImE|((0,fc) 


mico) 


(D4) 


S[ImE] is known as the Shannon-Jaynes entropy and the funetion m{co) is the so-ealled “default 
model”. In ease of no available data k) or infinitely large error a„ ^(M, k), the MEM 

proeedure will just give m{co) for ImE-[-(( 0 , k) beeause this funetion maximizes /’[ImE|//]. The 
default model ean thus be utilized to ineorporate already known information about ImE|(( 0 ,fc) 
into the analysis. 
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Collecting all the terms discussed above, we reach the final form of the probability /’[ImE|D//]: 

P[lmL\DH] oc P[D\lmLH]P[lmL\H] 

= ge[imq^ (D5) 

2[ImE] = aS[ImE] — L[ImE]. 

It is now merely a numerical problem to obtain the form of ImE|((0,fc) that maximizes 2[ImE] 
and is therefore the most probable ImE|((0, k) for given k) and H. For this task, we will 

use the Bryan algorithm [51]. 

Once ImE^ k) maximizing 2[ImE] for a fixed value of a is found, it is integrated out by 
averaging ImEj^ |((0, k) over a certain range of a, which then leads to our final result. Explicit 
formulae for this step and all other practical details specific to the application of MEM to QCD 
sum rules are discussed in [36, 44]. 

As a final point, let us mention here that Eqs.(27) and (28) give two independent sum rules, 
which have to be combined in the analysis of this work. How this can be done is explained in [44]. 
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